
# --- Libraries --- #
#devtools::install_github('jdwilson4/multilayer_extraction')
library(network)
library(MultilayerExtraction, quietly = TRUE)
library(Matrix, quietly = TRUE)
library(foreach, quietly = TRUE)  
library(doParallel, quietly = TRUE) 
library(parallel, quietly = TRUE)  
library(igraph, quietly = TRUE)
library(cccd)
library(countrycode)
library(xergm)
set.seed(1912)

setwd("YOUR WORKING DIRECTORY")


Mode <- function(x) {
  ux <- unique(x)
  ux[which.max(tabulate(match(x, ux)))]
}




# --- Multilayer Extraction --- #

# --- Weak signals --- #

vo_cands <- readRDS("data/comm_detection/candidate_layers/weak/vo_cands.rds") # always the same


dir <- "data/comm_detection/candidate_layers/weak/"
file_list <- list.files(dir, pattern = "sp_cands")


for(j in 1:length(file_list)){

  file.i <- file_list[[j]]
  
  print(j)
  
  sp_cands <- readRDS(paste(dir, file.i, sep = "")) # different embeddings



multi_comms_w_membs <- list()
multi_comms_w <- list()
for(i in 1970:1990){

  print(i)
  
  # M-E doesn't accomodate nodes that are isolates in each layer; so remove and add back at the end
  
  check_isos <- Reduce("+", list(sp_cands[[i]], vo_cands[[i]])) 
  check_isos2 <- check_isos[rowSums(check_isos)>0, colSums(check_isos)>0]
  
  s.adj <- adjust(sp_cands[[i]], check_isos2, remove = T)
  v.adj <- adjust(vo_cands[[i]], check_isos2, remove = T)
  
  node_names <- rownames(s.adj)
  
  rownames(s.adj) <- NULL; colnames(s.adj) <- NULL
  rownames(v.adj) <- NULL; colnames(v.adj) <- NULL
  
  l <- list(speeches = get.adjacency(graph_from_adjacency_matrix(s.adj)),
            votes = get.adjacency(graph_from_adjacency_matrix(v.adj)))
  
  n <- dim(l[[1]])[1]
  m <- length(l)
  
  net <- adjacency.to.edgelist(l)
  
  set.seed(3)
  comm_object <- multilayer.extraction(adjacency = net, seed = 3, min.score = 0, prop.sample = .10)
  obj <- refine(comm_object, k = Mode(comm_object$Diagnostics$Number.Communities), m = m, n = n)
  
  rownames(obj$Vertices) <- node_names
  m.v <- as.matrix(obj$Vertices)
  g.m <- m.v %*% t(m.v)
  
  g.m <- adjust(g.m, sp_cands[[i]], add=T, value=0)
  
  diag(g.m) <- 0
  
  multi_comms_w[[i]] <- g.m
  
  colnames(m.v) <- as.character(1:ncol(m.v))
  
  multi_comms_w_membs[[i]] <- m.v
  
}

saveRDS(multi_comms_w, file = paste("data/comm_detection/detected_comms/weak/weak_ties_",gsub("sp_cands_", "", file.i), sep = ""))
saveRDS(multi_comms_w_membs, file = paste("data/comm_detection/detected_comms/weak/weak_membs_",gsub("sp_cands_", "", file.i), sep = ""))


}











# --- Strong signals --- #

# --- Master node set:

sci.list <- readRDS("data/comm_detection/candidate_layers/strong/sci_masternodes.rds") # master node set for system-level analysis (i.e. any state who had a cooperation agreement)
mil.list <- readRDS("data/comm_detection/candidate_layers/strong/mil_masternodes.rds")
com.list <- readRDS("data/comm_detection/candidate_layers/strong/com_masternodes.rds")
fish.list <- readRDS("data/comm_detection/candidate_layers/strong/fish_masternodes.rds")
telco.list <- readRDS("data/comm_detection/candidate_layers/strong/telco_masternodes.rds")



prop <- c(.20, .25, .30) 

for(q in 1:length(prop)){

multi_comms_s_membs <- list()
multi_comms_s <- list()
for(i in 1970:1990){
  
  print(i)
  
  check_isos <- Reduce("+", list(sci.list[[i]], mil.list[[i]], com.list[[i]], fish.list[[i]], telco.list[[i]])) 
  check_isos2 <- check_isos[rowSums(check_isos)>0, colSums(check_isos)>0]
  
  sci.adj <- xergm.common::adjust(sci.list[[i]], check_isos2, remove = T)
  mil.adj <- xergm.common::adjust(mil.list[[i]], check_isos2, remove = T)
  com.adj <- xergm.common::adjust(com.list[[i]], check_isos2, remove = T)
  fish.adj <- xergm.common::adjust(fish.list[[i]], check_isos2, remove = T)
  telco.adj <- xergm.common::adjust(telco.list[[i]], check_isos2, remove = T)
  
  node_names <- rownames(sci.adj)
  
  rownames(sci.adj) <- NULL; colnames(sci.adj) <- NULL
  rownames(mil.adj) <- NULL; colnames(mil.adj) <- NULL
  rownames(com.adj) <- NULL; colnames(com.adj) <- NULL
  rownames(fish.adj) <- NULL; colnames(fish.adj) <- NULL
  rownames(telco.adj) <- NULL; colnames(telco.adj) <- NULL
  
  l <- list(sci = get.adjacency(graph_from_adjacency_matrix(sci.adj)),
            mil = get.adjacency(graph_from_adjacency_matrix(mil.adj)),
            com = get.adjacency(graph_from_adjacency_matrix(com.adj)),
            fish = get.adjacency(graph_from_adjacency_matrix(fish.adj)),
            telco = get.adjacency(graph_from_adjacency_matrix(telco.adj)))  
  
  n <- dim(l[[1]])[1]
  m <- length(l)
  
  net <- adjacency.to.edgelist(l)
  
  set.seed(3)
  comm_object <- multilayer.extraction(adjacency = net, seed = 3, min.score = 0, prop.sample = prop[q])
  obj <- refine(comm_object, k = Mode(comm_object$Diagnostics$Number.Communities), m = m, n = n)
  
  rownames(obj$Vertices) <- node_names
  
  m.v <- as.matrix(obj$Vertices)
  g.m <- m.v %*% t(m.v)
  
  g.m <- xergm.common::adjust(g.m, sci.list[[i]], add=T, value=0)
  
  diag(g.m) <- 0
  
  multi_comms_s[[i]] <- g.m
  
  colnames(m.v) <- as.character(1:ncol(m.v))
  
  multi_comms_s_membs[[i]] <- m.v
  
}

saveRDS(multi_comms_s, file = paste("data/comm_detection/detected_comms/strong/strong_ties_masternodes_prop",prop[q],".rds", sep = ""))
saveRDS(multi_comms_s_membs, file = paste("data/comm_detection/detected_comms/strong/strong_membs_masternodes_prop",prop[q],".rds", sep = ""))


}








# --- Weak node set:

sci.list <- readRDS("data/comm_detection/candidate_layers/strong/sci_weaknodes.rds") # weak node set (i.e. any state who spoke+voted in UN that year)
mil.list <- readRDS("data/comm_detection/candidate_layers/strong/mil_weaknodes.rds")
com.list <- readRDS("data/comm_detection/candidate_layers/strong/com_weaknodes.rds")
fish.list <- readRDS("data/comm_detection/candidate_layers/strong/fish_weaknodes.rds")
telco.list <- readRDS("data/comm_detection/candidate_layers/strong/telco_weaknodes.rds")


prop <- c(.20, .25, .30) 

for(q in 1:length(prop)){
  
  multi_comms_s_membs <- list()
  multi_comms_s <- list()
  for(i in 1970:1990){
    
    print(i)
    
    check_isos <- Reduce("+", list(sci.list[[i]], mil.list[[i]], com.list[[i]], fish.list[[i]], telco.list[[i]])) 
    check_isos2 <- check_isos[rowSums(check_isos)>0, colSums(check_isos)>0]
    
    sci.adj <- xergm.common::adjust(sci.list[[i]], check_isos2, remove = T)
    mil.adj <- xergm.common::adjust(mil.list[[i]], check_isos2, remove = T)
    com.adj <- xergm.common::adjust(com.list[[i]], check_isos2, remove = T)
    fish.adj <- xergm.common::adjust(fish.list[[i]], check_isos2, remove = T)
    telco.adj <- xergm.common::adjust(telco.list[[i]], check_isos2, remove = T)
    
    node_names <- rownames(sci.adj)
    
    rownames(sci.adj) <- NULL; colnames(sci.adj) <- NULL
    rownames(mil.adj) <- NULL; colnames(mil.adj) <- NULL
    rownames(com.adj) <- NULL; colnames(com.adj) <- NULL
    rownames(fish.adj) <- NULL; colnames(fish.adj) <- NULL
    rownames(telco.adj) <- NULL; colnames(telco.adj) <- NULL
    
    l <- list(sci = get.adjacency(graph_from_adjacency_matrix(sci.adj)),
              mil = get.adjacency(graph_from_adjacency_matrix(mil.adj)),
              com = get.adjacency(graph_from_adjacency_matrix(com.adj)),
              fish = get.adjacency(graph_from_adjacency_matrix(fish.adj)),
              telco = get.adjacency(graph_from_adjacency_matrix(telco.adj)))  
    
    n <- dim(l[[1]])[1]
    m <- length(l)
    
    net <- adjacency.to.edgelist(l)
    
    set.seed(3)
    comm_object <- multilayer.extraction(adjacency = net, seed = 3, min.score = 0, prop.sample = prop[q])
    obj <- refine(comm_object, k = Mode(comm_object$Diagnostics$Number.Communities), m = m, n = n)
    
    rownames(obj$Vertices) <- node_names
    
    m.v <- as.matrix(obj$Vertices)
    g.m <- m.v %*% t(m.v)
    
    g.m <- xergm.common::adjust(g.m, sci.list[[i]], add=T, value=0)
    
    diag(g.m) <- 0
    
    multi_comms_s[[i]] <- g.m
    
    colnames(m.v) <- as.character(1:ncol(m.v))
    
    multi_comms_s_membs[[i]] <- m.v
    
  }
  
  saveRDS(multi_comms_s, file = paste("data/comm_detection/detected_comms/strong/strong_ties_weaknodes_prop",prop[q],".rds", sep = ""))
  saveRDS(multi_comms_s_membs, file = paste("data/comm_detection/detected_comms/strong/strong_membs_weaknodes_prop",prop[q],".rds", sep = ""))
  
  
}





